ניתוח אזורים סטטיסטיים לאשכולות לפי גיל ומין
ד"ר עדו קליין
04/03/2021
ניתוח אשכולות לישובים ואזורים סטטיסטיים בישראל, לפי מין וגיל
מבוא
כולנו מניחים את ההנחות המוקדמות לגבי כל מגזר אוכלוסייה בישראל, במיוחד בימים אלו. אך במידה והיינו עוורים למגזר, האם ניתן לקבץ צורות יישוב שונות לפי אופי שונה של יישוב?
התשובה יכולה להינתן די בקלות על ידי ניתוח אשכולות באמצעות שפת התכנות הסטטיסטית R, על ידי ביצוע של הניתוח על נתוני גיל ומגדר של יישובים ואזורים סטטסיטיים.
להלן אציג את הניתוח המוצע על ידי.
שימו לב - הניתוח נהיה טכני קמעה לפעמים. זאת מפני שפוסט זה מציב לעצמו שתי מטרות - גם להציג את יכולתיה של שפת התכנות R לקהל דובר העברית, אך גם לבצע ניתוח נתונים מעניין. במידה ואתם משתעממים ממניפולציה של טבלאות, מציע לדלג על קטעים עם רקע בצבע הזה.
מתודולוגיה, או - תכנית עבודה
מה בתכנית?
1. טעינה וניקוי של הנתונים
2. ניתוח אשכולות
3. יצירת פירמידות גילאים
4. הצגת הנתונים על גבי מפה
שתי הערות טכניות לפני שמתחילים:
אני לא מאמין בהורדת המידע למחשב, אלא שימוש בקבצים זמניים בהורדה ישירה מהמקור. מה שאומר שאם כל הספריות הללו מותקנות לכם על המחשב, אתם יכולים להריץ הכל אצלכם ולראות שאינני מחרטט אתכם.
אני אסיר תודה ליונתן גת שהראה כאן איך ליישר עברית בrmarkdown.
טעינת הספריות הרלוונטיות
לניתוח זה נשתמש בחבילות הבאות:
{tidyverse}, המשמשת לעיבוד ו“שינוי צורה” של נתונים.
{httr}, המשמשת לעבודה עם כתובות אתרים ובקשות מידע מהם.
{readxl}, המשמשת לקריאת קבצי אקסל לתוך R.
{cluster}, המאפשרת ניתוח אשכולות במגוון כלים.
{sf}, המאפשרת קריאה ומניפולציה של קבצים המכילים נתונים מרחביים.
{mapview}, המאפשרת הצגה ויזאולית של נתונים מרחביים.
{leafpop}, המאפשרת הצגה של תמונות כפופאפים למפות.
{rmapshaper}, למניפולציה של ישויות מרחביות.
{randomcoloR}, לבחירת צבעים רנדומליים.
library(tidyverse)
library(httr)
library(readxl)
library(cluster)
library(sf)
library(mapview)
library(leafpop)
library(rmapshaper)
library(htmltools)
library(randomcoloR)טעינת הנתונים
אנו נסתמך על שני מאגרי נתונים - מאגר נתונים אחד אותו ננתח, ומאגר אחד שיסייע בהצגת הנתונים.
המאגר לניתוח הנתונים נמצא בדף אוכלוסייה ומרכיבי גידול ביישובים ובאזורים סטטיסטיים, 2019-2014, ואנו נשתמש בלוח 5 של שנת 2019.
ניתן לראות כי בקובץ של הלמ“ס מופיעות מספר עמודות רלוונטיות - זכרים, נקבות, ואיחודי א”ס 2019.
נוריד את הנתונים, ונביט בטבלת הזכרים:
url1 <- "https://www.cbs.gov.il/he/publications/doclib/2017/population_madaf/population_madaf_2019_7.xlsx"
GET(url1, write_disk(tf <- tempfile(fileext = ".xlsx")))## Response [https://www.cbs.gov.il/he/publications/doclib/2017/population_madaf/population_madaf_2019_7.xlsx]
## Date: 2021-03-04 21:54
## Status: 200
## Content-Type: application/vnd.openxmlformats-officedocument.spreadsheetml.sheet
## Size: 1.21 MB
## <ON DISK> C:\Users\User\AppData\Local\Temp\RtmpSayrHV\file347421d672dc.xlsxdf_male <- read_excel(tf, sheet =4,range = "A12:W2854")
df_female <- read_excel(tf, sheet =5,range = "A12:W2854")
df_unify <- read_excel(tf, sheet =6,range = "A9:D203")
df_male## # A tibble: 2,842 x 23
## `סמל צורת יישוב` `סמל יישוב` `שם יישוב` `אזור סטטיסטי` `אוכלוסייה בסוף~ גיל
## <chr> <dbl> <chr> <chr> <chr> <chr>
## 1 <NA> NA <NA> <NA> <NA> 4-0
## 2 "סה\"כ אוכלוסיי~ NA <NA> <NA> 4537918 4714~
## 3 "סה\"כ אוכלוסיי~ NA <NA> <NA> 4501377 4650~
## 4 "120" 3000 ירושלים סך הכל 467359 60301
## 5 "120" 3000 ירושלים 111 1444 242
## 6 "120" 3000 ירושלים 112 1442 177
## 7 "120" 3000 ירושלים 113 1992 444
## 8 "120" 3000 ירושלים 114 1952 433
## 9 "120" 3000 ירושלים 115 2319 546
## 10 "120" 3000 ירושלים 116 3539 824
## # ... with 2,832 more rows, and 17 more variables: ...7 <chr>, ...8 <chr>,
## # ...9 <chr>, ...10 <chr>, ...11 <chr>, ...12 <chr>, ...13 <chr>,
## # ...14 <chr>, ...15 <chr>, ...16 <chr>, ...17 <chr>, ...18 <chr>,
## # ...19 <chr>, ...20 <chr>, ...21 <chr>, ...22 <chr>, ...23 <chr>ניתן לראות כי למרות שביקשנו טווח מוגדר (R לא יודע לנחש את מיקום הטבלה כל פעם), חלק משמות העמודות היו בעייתים. לכן יש צורך בסידור הנתונים (ועל הדרך אנחנו גם ננקה אותם).
סידור וניקוי
כפי שניתן לראות, R לא אוהב שמות עמודות בעברית או שמות עמודות כפולות בתאים ממוזגים, ולכן יש צורך בשינוי שמות העמודות.
בנוסף, זרקנו כמה עמודות זבל מהתהליך ועמודות שלא קשורות ליישובים - 3 הראשונות.
כל פעם שהמילה סך הכל הופיעה במספר האזור הסטטיסטי עבור ערים - השמטנו את התצפית על מנת להימנע מכפילות.
ניקינו את הנתונים כך שכל התאים עם ערכי האוכלוסיה יהיו מספריים, ושתאים עם ערך לא ידוע בעקבות ההמרה יהיו שווים 0.
לבסוף, איחדנו בין טבלת הזכרים והנקבות, סינון אזורים בהם סכום הגברים והנשים הוא 0, ומחקנו מלל מהעמודה המצביעה על שם האזור.
# setting workable column name
alt_col_names <- c("type_sym","setl_sym","setl_name","tract","pop2019")
ages <- as.list(as.data.frame(t(df_male[1,])))[[1]]
ages <- ages[!is.na(ages)]
ages <- str_replace(ages,"-","_")
ages <- str_replace(ages,"\\+","_")
ages <- paste0("Age_",ages)
# manual fix
ages[18] <- "Age_up_85"
alt_col_names <- c(alt_col_names,ages)
colnames(df_male) <- alt_col_names
colnames(df_female) <- alt_col_names
# discarding useless rows
df_male <- df_male[-c(1:3),]
df_female <- df_female[-c(1:3),]
# filtering sums
df_male <- df_male %>%
group_by(setl_sym) %>%
mutate(n = n()) %>%
ungroup() %>%
filter(tract != "סך הכל"|n == 1) %>%
select(-n)
df_female <- df_female %>%
group_by(setl_sym) %>%
mutate(n = n()) %>%
ungroup() %>%
filter(tract != "סך הכל"|n == 1) %>%
select(-n)
# convert relevant columns to numeric and fill NA with zero
df_male <- df_male %>%
mutate_at(vars(pop2019:Age_up_85),as.numeric) %>%
mutate_at(vars(pop2019:Age_up_85),replace_na,replace = 0)
df_female <- df_female %>%
mutate_at(vars(pop2019:Age_up_85),as.numeric) %>%
mutate_at(vars(pop2019:Age_up_85),replace_na,replace = 0)
# join males and females, filter zero pop, tract text to NA
df <- df_male %>%
left_join(df_female, by = c("type_sym", "setl_sym" ,"setl_name" ,"tract" ),suffix = c(".male",".female")) %>%
filter(pop2019.male + pop2019.female != 0) %>%
mutate(tract = ifelse(tract == "סך הכל", "1",tract))
df## # A tibble: 2,753 x 42
## type_sym setl_sym setl_name tract pop2019.male Age_4_0.male Age_9_5.male
## <chr> <dbl> <chr> <chr> <dbl> <dbl> <dbl>
## 1 120 3000 ירושלים 111 1444 242 205
## 2 120 3000 ירושלים 112 1442 177 133
## 3 120 3000 ירושלים 113 1992 444 352
## 4 120 3000 ירושלים 114 1952 433 351
## 5 120 3000 ירושלים 115 2319 546 472
## 6 120 3000 ירושלים 116 3539 824 648
## 7 120 3000 ירושלים 121 2398 240 199
## 8 120 3000 ירושלים 122 1580 196 137
## 9 120 3000 ירושלים 123 2311 261 188
## 10 120 3000 ירושלים 124 1098 97 101
## # ... with 2,743 more rows, and 35 more variables: Age_14_10.male <dbl>,
## # Age_19_15.male <dbl>, Age_24_20.male <dbl>, Age_29_25.male <dbl>,
## # Age_34_30.male <dbl>, Age_39_35.male <dbl>, Age_44_40.male <dbl>,
## # Age_49_45.male <dbl>, Age_54_50.male <dbl>, Age_59_55.male <dbl>,
## # Age_64_60.male <dbl>, Age_69_65.male <dbl>, Age_74_70.male <dbl>,
## # Age_79_75.male <dbl>, Age_84_80.male <dbl>, Age_up_85.male <dbl>,
## # pop2019.female <dbl>, Age_4_0.female <dbl>, Age_9_5.female <dbl>,
## # Age_14_10.female <dbl>, Age_19_15.female <dbl>, Age_24_20.female <dbl>,
## # Age_29_25.female <dbl>, Age_34_30.female <dbl>, Age_39_35.female <dbl>,
## # Age_44_40.female <dbl>, Age_49_45.female <dbl>, Age_54_50.female <dbl>,
## # Age_59_55.female <dbl>, Age_64_60.female <dbl>, Age_69_65.female <dbl>,
## # Age_74_70.female <dbl>, Age_79_75.female <dbl>, Age_84_80.female <dbl>,
## # Age_up_85.female <dbl>עיבוד מקדים לנתונים
על מנת לייצר בסיס נתונים להזין אותו לפונקציה שמבצעת את אלגוריתם האשכולות, יש לנרמל את הנתונים עבור כל שכבת גיל, לבטא אותם בצורה יחסית לכלל האוכלוסיה באזור המסוים, ולשמור רק את עמודות הגיל.
df_clus <- df %>%
mutate(pop2019 = pop2019.male + pop2019.female) %>%
mutate_at(vars(starts_with("Age")),~./pop2019) %>%
select(starts_with("Age"))
df_clus## # A tibble: 2,753 x 36
## Age_4_0.male Age_9_5.male Age_14_10.male Age_19_15.male Age_24_20.male
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.0817 0.0692 0.0729 0.0405 0.0419
## 2 0.0615 0.0462 0.0684 0.0531 0.0698
## 3 0.113 0.0899 0.0536 0.0184 0.0232
## 4 0.110 0.0895 0.0553 0.0145 0.0196
## 5 0.118 0.102 0.0515 0.0149 0.0140
## 6 0.117 0.0923 0.0433 0.0135 0.0197
## 7 0.0483 0.0400 0.0436 0.0378 0.0332
## 8 0.0605 0.0423 0.0364 0.0290 0.0318
## 9 0.0540 0.0389 0.0300 0.0259 0.0300
## 10 0.0425 0.0442 0.0469 0.0333 0.0285
## # ... with 2,743 more rows, and 31 more variables: Age_29_25.male <dbl>,
## # Age_34_30.male <dbl>, Age_39_35.male <dbl>, Age_44_40.male <dbl>,
## # Age_49_45.male <dbl>, Age_54_50.male <dbl>, Age_59_55.male <dbl>,
## # Age_64_60.male <dbl>, Age_69_65.male <dbl>, Age_74_70.male <dbl>,
## # Age_79_75.male <dbl>, Age_84_80.male <dbl>, Age_up_85.male <dbl>,
## # Age_4_0.female <dbl>, Age_9_5.female <dbl>, Age_14_10.female <dbl>,
## # Age_19_15.female <dbl>, Age_24_20.female <dbl>, Age_29_25.female <dbl>,
## # Age_34_30.female <dbl>, Age_39_35.female <dbl>, Age_44_40.female <dbl>,
## # Age_49_45.female <dbl>, Age_54_50.female <dbl>, Age_59_55.female <dbl>,
## # Age_64_60.female <dbl>, Age_69_65.female <dbl>, Age_74_70.female <dbl>,
## # Age_79_75.female <dbl>, Age_84_80.female <dbl>, Age_up_85.female <dbl>הרצת אלגוריתם לניתוח אשכולות
כעת, כאשר הנתונים מסודרים כיאות, ניתן סוף סוף להריץ את הניתוח הרצוי.
קיימים מספר אלגוריתמים לביצוע הניתוח, אנו נשתמש בניתוח אשכלות היררכי על בסיס חלוקה, היינו - בשלב הראשון כל התצפיות מקובצות לאשכול אחד, ולאט לאט הן מתפרקות לאשכולות שונים, על בסיס המרחק שלהן אחת מהשנייה - ככל שתצפית רחוקה יותר מחברותיה, כך היא תתפרק מהאשכול מוקדם יותר. לבסוף, מקבלים חלוקה של כל התצפיות, לפי סוגים של אזורים המתחברים אחד לשני על סמך מרחק. עוד על שיטת החלוקה - בנספח.
z <- diana(df_clus)
plot(z, main="plotit", which.plot = 2)
k <- 44
clus_table <- table(cutree(z,k = k)) %>% data.frame()
clusters <- cutree(z,k = 44)הניתוח שביצעתי מראה שכאשר מחלקים אותו ל44 חלקים, מקבלים אשכולות “מעניינים”. הגרף הנ"ל נקרא דנדוגרמה, והוא נועד להראות מערכות יחסים היררכיות בין אובייקטים. במקרה שלנו, הוא מראה אינטואיטיבית באיזה מרחק האשכולות נפרדים אחד מהשני.
לא ארחיב כאן בהסבר הטכני, ואשלח את הקוראים לבדוק עוד כאן.
השבת האשכולות אל הטבלה
הפלט של התהליך יוצר לנו וקטור המשייך כל תצפית לאשכול מסוים. נחבר את האשכול חזרה אל הטבלה המקורית, ונבחן את התפלגות האשכולות לפי סוג.
df$clusters <- clusters
clus_table %>% ggplot(aes(Var1,Freq)) + geom_col() + labs(x="מספר אשכול",y="כמות")
df <- df %>%
group_by(clusters) %>%
mutate(n = n()) %>%
ungroup() %>%
mutate(clusters = ifelse(n < 20, 99,clusters)) %>%
select(-n)ניתן לראות שאנו מקבלים אשכול אחד שלוקח את רוב הדגימות, אשכול אחד בגודל בינוני, עוד ארבע אשכולות קטנים והשאר זעירים. לכן, בהמשך הניתוח, נתייחס לאשכולות הזעירים כאל אחד, שכן הם מבטא מקרים אקזוטיים, וככל הנראה מדובר ביישובים קטנים מאוד/מוסדות.
יצירת פירמידת גילאים
לאחר שביצענו את הניתוח המסווג את הישובים והאזורים הסטטיסטים לאשכולות, נבצע ויזואליזציה של פירמידות הגילאים באזורים השונים.
הויזואליזציה תתחלק לשניים: ויזואליזציה של האשכולות המייצגים, ויזואליזציה עבור כל אזור.
כאן יתבצע גם תהליך של שיום של האשכולות.
הקוד לשרטוט פירמידת הגילאים מבוסס על התשובה הזאת.
eamp <- df %>%
nest(starts_with("Age")) %>%
mutate(text = paste(setl_name,"אזור:",tract),
plo = map2(data,text,~t(.x) %>%
as.data.frame() %>%
magrittr::extract(1) %>%
rownames_to_column("age_split") %>%
separate(age_split,c("del","to","from_dirty"),sep = "_") %>%
separate(from_dirty,c("from","gender"),sep = "\\.") %>%
unite("age_group",from,to, sep = "-",remove = FALSE) %>%
mutate(levels = factor(age_group,levels(fct_reorder(unique(age_group),1:length(unique(age_group)))))) %>%
select(-del) %>%
ggplot(aes(x= ifelse(gender == "male", -V1,+V1),y = levels, fill = gender))+
geom_col() +
lemon::scale_x_symmetric(labels = abs) +
labs(x = "Population",y = "age group",title = .y) +
theme(plot.title = element_text(hjust = 0.5))
))
facet_plot <- df %>%
group_by(clusters) %>%
mutate(n = n()) %>%
group_by(clusters,n) %>%
summarise_at(vars(pop2019.male,pop2019.female, starts_with("Age")), sum) %>%
mutate(pop2019 = pop2019.male + pop2019.female) %>%
select(-pop2019.male,-pop2019.female) %>%
mutate_at(vars(starts_with("Age")),~./pop2019) %>%
gather(key,value,-clusters,-pop2019,-n) %>%
separate(key,c("del","to","from_dirty"),sep = "_") %>%
separate(from_dirty,c("from","gender"),sep = "\\.") %>%
unite("age_group",from,to, sep = "-") %>%
select(-del) %>%
mutate(levels = factor(age_group,levels(fct_reorder(unique(age_group),1:length(unique(age_group)))))) %>%
filter(clusters != 99) %>%
ggplot(aes(x= ifelse(gender == "male", -value,+value),y = levels, fill = gender))+
geom_col() +
geom_label(aes(x= -0.05 , y = 14, label = paste("n:",n)),fill = "white") +
geom_label(aes(x= -0.075 , y = 17, label = paste("tot_pop:",pop2019)),fill = "white",size = 4) +
lemon::scale_x_symmetric(labels = abs) +
labs(x = "Fraction of population",y = "age group") +
facet_wrap(~clusters)
facet_plot 
כפי שניתן לראות, התפלגות הגילאים באשכולות השונים שונה בכל אחד מהאשכולות. ניתן לתת שם לכל אחד מהאשכולות השונים:
אשכול מספר 1 מזכיר מאוד פירדמידות גילאים המאפיינות מדינות עולם מתפתח, ולכן נקרא לאשכול זה - “עולם מתפתח”.
אשכול מספר 2 בולט בכך שהוא מכיל שכבה גדולה של הורים צעירים לילדים רבים. לכן נקרא לו “משפחות צעירות ברוכות ילדים”.
אשכול מספר 3 מזכיר מאוד פירמידות גילאים המאפיינות מדינות עולם מפותח, לכן נקרא לו “מאוזן”.
אשכול מספר 4 מכיל משפחות מעט יותר מבוגרות מאשכול מספר 2 המכילות יחסית פחות ילדים, לכן נקרא לו “משפחות בינוניות”.
אשכול מספר 5 מכיל משפחות מבוגרות יותר מאשכול 4, לכן נקרא לו “משפחות מבוגרות”.
אשכול מספר 6 הוא אחד ההיילייטים של הניתוח. לאחר שגיליתי את קיומו של אשכול זה, הייתי צריך לחשוב מה אני רואה. לאחר ששלחתי מייל ללמ"ס כדי להבהיר מספר דברים, נאמר לי כי:
תשובה זו אוששה את חשדי לגבי מהות אשכול זה. מדובר באשכול המלא בפנימיות לבנים, תוך כדי הימצאות של פירמידת גילאים של עולם מתפתח. לכן לאשכול זה נקרא “ישובים עם ישיבה”.
אשכול מספר 10 בולט בכך שהוא מכיל אוכלוסיה צעירה מרובה ויחסית מעט ילדים, לכן נקרא לו “צעירים”.
אשכול 22 הוא מקרה קיצוני של אשכול 11, לכן נקרא לו “מלא צעירים”. שימו לב שלכל אחד מהאשכולות מצורפים כמות האשכולות ומספר האנשים הכולל המתגורר באשכול.
clusters_names <- c("עולם מתפתח",
"משפחות צעירות ברוכות ילדים",
"מאוזן",
"משפחות בינוניות",
"משפחות מבוגרות",
"ישובים עם ישיבה",
rep(NA,3),
"צעירים",
rep(NA,11),
"ים צעירים")
facet_plot$data <- facet_plot$data %>%
mutate(clusters = clusters_names[clusters])
facet_plot
הורדת מידע גיאוגרפי
כעת, לאחר שסיימנו לנתח את האשכולות השונים, חשבתי שיהיה מעניין לראות אותם על מפה, כדי להבחין בדפוסים מרחביים כאלה ואחרים. בשביל לעשות את זה, נוריד את שכבת האזורים הסטטיסטיים והישובים של הלמ"ס מכאן.
ניקוי מידע גיאוגרפי, וחיבורו ליתר הנתונים
כאשר הורדתי את המידע הגיאוגרפי גיליתי כי הפילוח לפי גילאים היה בתוך השכבה, אם כי לא היה פילוח לפי מגדר. לכן, במידה מסויימת, יתכן שחלק מסידור הנתונים היה מיותר. אם כי, במידה ולא היינו מנתחים את האקסל, לא היינו מוצאים את אשכול מספר 6….
כשניסיתי להנגיש את השכבה בצורה אינטראקטיבית, גיליתי שהיא כבדה מדי, ולכן העברתי אותה תהליך של פישוט.
שוב, מאוד טכני - פישוט שמתבצע על הרבה ישויות שחולקות גבול משותף הוא בעייתי, כי אם הוא מתבצע בנפרד הוא עלול לגרום לחפיפות ולחללים בין הישויות. לכן השתמשתי בחבילה המאפשרת פישוט תוך כדי שמירה על הגבולות בין הישויות. עוד על כך בלינק הקודם, חנונים!
כעת נתחיל בניקוי של השכבה הגאוגרפית. איזורים בהם ערך האוכלוסיה הוא NA, יסוננו החוצה מהשכבה.
בנוסף, נבצע חיבור טבלאי של השכבה לטבלאות הקודמות שיצרנו - נתונים וגרפים של פירמידת הגילאים לכל אזור/ישוב.
shp <- st_cast(shp,"MULTIPOLYGON")
shp <- ms_simplify(shp,keep_shapes = TRUE)
maplot <- shp %>%
filter(!is.na(Pop_Total),!is.na(STAT11)) %>%
mutate(STAT11 = as.character(STAT11)) %>%
left_join(df, by =c("SEMEL_YISHUV" = "setl_sym", "STAT11" = "tract")) %>%
left_join(eamp %>% select(setl_sym,tract,plo), by =c("SEMEL_YISHUV" = "setl_sym", "STAT11" = "tract")) %>%
filter(!is.na(clusters)) %>%
mutate(clusters = clusters_names[clusters],
clusters = ifelse(is.na(clusters), "אקזוטי",clusters),
plo = map2(plo,clusters,~.x + labs(subtitle = .y) + theme(plot.subtitle = element_text(hjust = 0.5))))הצגת הנתונים על גבי מפה
ועכשיו לרגע שכולנו חיכינו לו - המפה! המפה מראה את גבולות האיזורים ומסמנת לפי צבע את האשכול אליו הם שייכים. בנוסף, לחיצה על איזור מאפשרת לראות את פירמידת הגילאים שלו. הערה קטנה לגבי התנחלויות: למעט הערים בהתנחלויות, הלמ"ס לא שומר על שכבות של הגבולות של כל אחד מהישובים, אלא מסמן אותם כנקודה שעליה עושים חיץ קטן. אם חשקה נפשכם לראות נתונים של התנחלות כזו או אחרת, פשוט תתמרכזו אליה.
מזמין אתכם לגלות דברים מעניינים. למשל - האם אשכולות מסוימים נוטים להיות קרובים אחד לשני? האם אתם יכולים להזדהות עם הניתוח לגבי אזורים מסוימים בעיר שלכם? האם המפה מחדשת לכם דברים שלא ידעתם?
בכל מקרה, חובה עליכם לעשות זום אין על אזור תל אביב. מאוד מעניין מה שקורה שם.
mapviewOptions(vector.palette = randomColor(9))
map <- mapview(maplot,zcol = "clusters",popup = popupGraph(maplot$plo))
map@map